#coding=utf-8
'''
Created on 2011-3-14

@author: zarra

基于差分的R波检测
'''
import mitEcgDbReader
import matplotlib as mpl
mpl.rcParams['font.sans-serif'] = ['SimHei']
mpl.rcParams['axes.unicode_minus'] = False
from pylab import *
import numpy as np
import scipy.signal as signal




ecgData = mitEcgDbReader.EcgDataFile('../data/100.hea')

I=0.2 #阈值
S= float(ecgData.head.sampleFreq) #采样律
D= S * 0.2 #规定的最小 R-R间距
    
channel = ecgData.readAll()[0]
x=channel.real_data
n=len(x)
xaxis = np.linspace(0,n/S,n)
#x(n-2)
x_2= np.linspace(0, 0, n)
x_2[:-2]=x[2:]

#x(n-4)
x_4= np.linspace(0, 0, n)
x_4[:-4]=x[4:] 

#y0=|x(n)-x(n-2)|
y0=np.abs(x-x_2)
#y1=|x(n)-2x(n-2)+x(n-4)|
y1=np.abs(x-2*x_2+x_4)
#y2=1.3y0+1.1y1
y2=1.3*y0+1.1*y1

#y3为y2的8点平滑
#外加阈值判定
y3=(y2[:-8]+y2[1:-7]+y2[2:-6]+y2[3:-5]+\
    y2[4:-4]+y2[5:-3]+y2[6:-2]+y2[7:-1]+y2[8:])/8
y3_binarization = np.frompyfunc( lambda  value:  value if value > I else 0,1,1)
y3 = y3_binarization(y3)
y3_len = len(y3)
r_wave_begin = np.linspace(0,0,y3_len/2)
i=0
index = 0
while i<y3_len:
    if y3[i]>0:
        r_wave_begin[index ]=float(i)/S
        index+=1
        i+=D
        continue
    i+=1

r_wave_begin = r_wave_begin[:index-1]
T= len(r_wave_begin)
STR1 = '心跳总数:%d'%T
r_wave_dst = r_wave_begin[1:]-r_wave_begin[:-1]
STR2 =  '平均R-R间隔:%f'%(float(sum(r_wave_dst))/len(r_wave_dst))

Tmax = 4000
figure(figsize=(8,4))
subplot(311)
plot(xaxis[:Tmax],x[:Tmax])
title(u"原心电波形 "+STR1+' '+STR2)

subplot(312)
plot(xaxis[:Tmax],y2[:Tmax])
title(u"y2波形")

subplot(313)

plot(xaxis[:Tmax],y3[:Tmax])
title(u"y3波形")
show()
if __name__ == '__name__':
    pass

    
    



